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Abstract  We  conduct  discrete  spectrum  analyses  for  a 
selection  of  mixed  discretization  schemes  for  the  Stokes 
eigenproblem.  In  particular,  we  consider  the  MINI  el¬ 
ement,  the  Crouzeix-Raviart  element,  the  Marker-and- 
Cell  scheme,  the  Taylor-Hoocl  element,  the  Qk/Pk-i  el¬ 
ement,  the  divergence-conforming  discontinuous  Galerkin 
method,  and  divergence-conforming  B-splines.  For  each 
of  these  schemes,  we  compare  the  spectrum  for  the  con¬ 
tinuous  Stokes  problem  with  the  spectrum  for  the  dis¬ 
crete  Stokes  problem,  and  we  discuss  the  relationship 
of  eigenvalue  errors  with  solution  errors  associated  with 
unsteady  viscous  flow  problems. 

Keywords  Stokes  eigenproblem  •  mixed  methods  • 
discrete  spectrum  analysis 

1  Introduction 

One  of  the  most  important  roles  of  numerical  analysis  is 
to  quantify  the  induced  error  of  a  particular  discretiza¬ 
tion  scheme.  Classical  a  priori  error  estimation  tech¬ 
niques  allow  one  to  identify  the  asymptotic  convergence 
rate  of  a  scheme,  but  they  rarely  give  much  insight 
into  a  method’s  pre-asymptotic  characteristics.  Alter¬ 
natively,  one  may  utilize  discrete  spectrum  analysis  to 
analyze  the  approximation  properties  of  all  the  scales  of 
a  discretization  scheme.  In  discrete  spectrum  analysis, 
one  directly  compares  the  spectrum  of  a  chosen  dif¬ 
ferential  operator  with  the  spectrum  of  the  discretized 
version  of  the  differential  operator.  Discrete  spectrum 
analysis  has  a  long  history  in  structural  dynamics  where 
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the  eigenvalues  of  a  linear  structural  system  correspond 
to  the  squares  of  the  natural  frequencies  [19].  For  lin¬ 
ear  structural  dynamics  problems,  it  can  be  shown  that 
errors  in  natural  frequency  lead  to  a  linear  accumu¬ 
lation  of  solution  error  in  time  which  scales  with  fre¬ 
quency  error.  For  nonlinear  structural  dynamics  prob¬ 
lems,  induced  coupling  of  low  and  high  modes  exacer¬ 
bates  this  accumulation  of  error,  ultimately  resulting  in 
non-robust  numerical  schemes.  Hence,  a  discretization 
technique  whose  discrete  spectrum  closely  matches  the 
exact  spectrum  is  highly  desirable.  Discrete  spectrum 
analysis  can  also  be  utilized  to  analyze  diffusive  prob¬ 
lems  where  the  eigenvalues  of  the  underly  differential 
systems  correspond  to  the  dissipation  rates  of  distinct 
eigenmodes  [2] .  While  it  can  be  easily  shown  that  errors 
in  dissipation  rates  do  result  in  enlarged  solution  error, 
these  errors  typically  decay  exponentially  in  time.  Nev¬ 
ertheless,  quantification  of  such  errors  is  important  for 
diffusive  problems  exhibiting  a  wide  array  of  spatial  and 
temporal  scales  such  as  viscous  turbulent  fluid  flow. 

In  this  paper,  we  conduct  a  discrete  spectrum  anal¬ 
ysis  for  a  selection  of  mixed  discretization  schemes  for 
the  Stokes  eigenproblem.  For  simplicity,  we  restrict  our 
attention  to  the  Stokes  eigenproblem  posed  on  the  two- 
dimensional  torus.  The  eigenvalues  of  the  Stokes  eigen¬ 
problem  correspond  to  the  dissipation  rates  of  the  natu¬ 
ral  modes  of  unsteady  Stokes  flow  as  well  as  the  squares 
of  the  natural  frequencies  of  linear  incompressible  elas- 
todynamics.  Furthermore,  in  the  context  of  incompress¬ 
ible  turbulent  fluid  flow,  these  eigenvalues  approximately 
match  the  dissipation  rates  of  modes  with  sufficiently 
high  wave  number,  namely  modes  in  the  viscous  range, 
and  comprise  a  significant  fraction  of  the  dissipation 
rates  of  modes  at  the  tail  end  of  the  inertial  range. 
Consequently,  a  scheme  with  favorable  resolution  prop- 
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erties  should  be  characterized  by  a  discrete  spectrum 
which  closely  mimics  the  exact  spectrum. 

An  outline  of  this  paper  is  as  follows.  In  Section 
2,  we  consider  the  unsteady  Stokes  flow  problem  and 
express  its  solution  in  terms  of  the  eigenvalues  and 
eigenmodes  of  the  Stokes  eigenproblem.  We  also  dis¬ 
cuss  numerical  approximation  of  the  unsteady  Stokes 
flow  problem  and  its  relationship  to  the  Stokes  eigen¬ 
problem,  and  we  discuss  estimation  of  solution  errors. 
In  Section  3,  we  conduct  discrete  spectrum  analyses  for 
a  variety  of  mixed  discretization  schemes  for  the  Stokes 
eigenproblem,  starting  with  the  MINI  element  in  Sub¬ 
section  3.1  and  ending  with  divergence-conforming  B- 
splines  in  Subsection  3.7.  In  Section  4,  we  draw  conclu¬ 
sions. 


2  Unsteady  Stokes  Flow  and  the  Stokes 
Eigenproblem 

To  properly  motivate  discrete  spectrum  analysis  and 
its  role  in  characterizing  errors  associated  with  time- 
dependent  problems,  we  consider  the  unsteady  Stokes 
flow  problem  in  this  section  and  express  its  solution 
in  terms  of  eigenvalues  and  eigenmodes.  For  simplicity, 
let  us  consider  the  unsteady  Stokes  problem  on  the  two- 
dimensional  torus  subject  to  no  applied  forcing.  To  this 
effect,  we  set  17  =  (0,27r)2  and  introduce  the  Sobolev 
spaces 

H  lerW)  :={u£  H1(17)  : 

u(-,0)  =  u(-,2tt),u(0,  •)  =  u(2tt,-)} 


and 


L2per(f2)  :=  L2{Q). 


lows. 


Find  u  G  VT,Per  and  p  G  Qr,Per  such  that  u(0)  = 
u0  and,  for  almost  every  t  G  (0,  T), 

(<9tu(f),v)  +  fc(u(i),v)  -  b(p(t),v) 

+b(q,u(t))  =  0 


(P){ 


for  all  v  G  Hper.(l7)  and  q  G  Lper(f2)  where 


fc(u,  v)  =  v  /  Vu  :  Vvdx, 
J  n 

b(p ,  v)  =  /  pV  •  vdx. 

J  n 


Above,  u  denotes  the  flow  velocity  of  a  fluid  moving 
through  the  domain  17,  p  denotes  the  pressure  acting  on 
the  fluid  divided  by  the  fluid  density,  and  v  denotes  the 
kinematic  viscosity  of  the  fluid.  Problem  (P)  consists  of 
an  evolutionary  equation  subject  only  to  dissipation  in 
time.  Now,  defining  the  space  of  divergence-free  func¬ 
tions  as 

H ler(n)  :=  {v  G  Hper(l7)  :  V  •  v  =  0}  , 


let  us  consider  the  following  eigenproblem. 


(  Find  u  G  Hper(l7)  and  A  G  R+  such  that 
=  1  and 


llullL2(i?) 

H  u,v)  =  A(u,v)L2 


l  2(n) 


{  for  all  v  G  H  (17). 


The  orthonormal  eigenmodes  corresponding  to  the  above 
problem  are  explicitly  known  and  comprise  the  set 


{2  = 


curl<?!> 


|  curb/)  1 1 


l  2(n) 


:  (j>  =  exp  (ik\X  +  ik^y) , 
hi,  fc2  G  Z,  k\  +  k'l  ^  0}  . 


Furthermore,  for  a  given  end-time  T  >  0,  we  define  the  We  enumerate  these  eigenmodes  and  their  correspond- 
spaces  ing  eigenvalues  as  {u„,An}^x  where 


VT,Per  :={vG  i2(0,  T;  Hper(17))  : 

dt.v  G  L2(0,T;  H*er(l7))} 

and 

Qr,Per  :=  T2(0,  T;  Lper(f2)) 

where  H*er(l7)  denotes  the  dual  space  of  Hper(l7).  Given 
u0  G  Hper(l7),  the  problem  of  interest  then  reads  as  fol- 


Ai  <  A2  <  •  •  • . 

Finally,  a  direct  calculation  shows  that  the  exact  veloc¬ 
ity  solution  of  (P)  can  be  written  as 

OO 

u(x,  t)  =  A0  +  ^2  A„u„(x)  exp(— A„t) 


where 


Uo(x)dx 
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Fig.  1  Degrees  of  freedom:  (a)  the  MINI  element,  (b)  the 
non-conforming  Crouzeix-Raviart  element,  (c)  the  Marker- 
and-Cell  finite  difference  scheme.  Velocity  DOF  are  on  the 
left  and  pressure  DOF  are  on  the  right. 


and 

An=  u0(x)  •  un(x)dx. 

J  n 

Hence,  the  velocity  solution  of  (P)  depends  only  on 
the  initial  condition  u0,  the  eigenmodes  {un}^Li,  and 
the  eigenvalues  {Anj^L-p  Furthermore,  the  eigenvalues 
gives  the  dissipation  rates  of  the  individual 
eigenmodes.  We  may  repeat  all  of  the  above  calcula¬ 
tions  at  the  discrete  level  to  express  our  discrete  ve¬ 
locity  solution  entirely  in  terms  of  the  discrete  initial 

r  ,}N 

condition  Uq,^,  discrete  eigenmodes  <  u„  >  ,  and  dis- 

l  J  n—1 

Crete  eigenvalues  {Ajj}  x  where  N  denotes  the  total 
number  of  discretely  divergence-free  modes  correspond¬ 
ing  to  nonzero  eigenvalues.  A  natural  question  to  then 
ask  is  how  well  do  the  discrete  eigenvalues  approximate 
the  first  N  exact  eigenvalues.  This  question  can  be  an¬ 


swered  to  some  degree  using  tools  arising  from  func¬ 
tional  analysis.  For  example,  it  is  well-known  that  the 
discrete  eigenvalues  arising  from  a  conforming,  stable, 
and  balanced  mixed  finite  element  approximation  of  the 
Stokes  eigenproblem  satisfy  the  relationship 

\\n-\hn\<Ch2k 

provided  the  product  nh  is  sufficiently  small  where  h  is 
the  mesh  size,  k  is  the  polynomial  degree  of  the  discrete 
velocity  space,  and  C  is  a  positive  constant  independent 
of  h  [4].  Unfortunately,  this  relationship  only  gives  an 
error  bound  for  discrete  eigenvalues  corresponding  to 
the  lowest  discrete  eigenmodes.  For  simple  eigenprob- 
lerns  such  as  those  arising  from  one-dimensional  struc¬ 
tural  dynamics,  one  can  employ  discrete  Fourier  anal¬ 
ysis  to  exactly  characterize  discrete  eigenvalue  errors 
across  the  entire  spectrum  [21],  though  this  analysis  is 
intractable  in  the  context  of  multi-dimensional  Stokes 
flow.  Hence,  in  this  paper,  we  resort  to  direct  computa¬ 
tion  to  compute  discrete  eigenvalue  errors  and  compare 
discrete  spectra  with  the  spectrum  corresponding  to  the 
continuous  Stokes  eigenproblem. 


3  Discrete  Spectrum  Analyses 

3.1  The  MINI  Element 

We  begin  our  discrete  spectrum  analyses  with  the  so- 
called  MINI  element.  Originally  introduced  by  Arnold, 
Brezzi,  and  Fortin  in  [1],  the  MINI  element  is  arguably 
the  most  efficient  mixed  discretization  scheme  for  the 
approximation  of  Stokes  flow.  The  velocity  space  of  the 
MINI  element  is  defined  to  be  that  of  piecewise  continu¬ 
ous  linears  enriched  by  the  space  of  element-wise  cubic 
bubbles  over  a  given  triangulation,  and  the  pressure 
space  is  defined  to  be  only  that  of  piecewise  continuous 
linears.  The  degrees  of  freedom  (DOF)  for  this  element 
are  illustrated  in  Fig.  1(a).  It  should  be  noted  that  if 
one  statically  condenses  the  cubic  bubbles  over  each 
element,  one  arrives  at  a  formulation  that  is  nearly  al¬ 
gebraically  equivalent  to  the  pressure-stabilized  Pi /Pi 
velocity /pressure  pair  [20]. 

We  have  computed  the  discrete  Stokes  spectrum 
corresponding  to  the  MINI  element  on  a  structured  tri¬ 
angular  mesh  with  1,568  elements  and  plotted  a  com¬ 
parison  of  the  discrete  spectrum  with  the  exact  spec¬ 
trum  in  Fig.  2.  In  the  figure,  n  denotes  the  given  mode 
number  while  N  denotes  the  total  number  of  discrete 
modes.  For  this  calculation,  N  =  3,  919.  Note  that  while 
the  initial  part  of  the  spectrum  is  fairly  well-resolved, 
the  discrete  eigenvalues  corresponding  to  n/N  >  0.15 
exhibit  more  than  fifty  percent  error.  As  the  discrete 
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Fig.  2  Stokes  spectrum  for  the  MINI  element. 


eigenvalues  are  larger  in  magnitude  than  the  continu¬ 
ous  eigenvalues,  this  indicates  that  in  an  unsteady  flow 
simulation,  eighty  five  percent  of  the  discrete  modes 
will  dissipate  much  faster  than  their  continuous  coun¬ 
terparts.  It  should  be  indicated  that  the  upper  part  of 
the  discrete  spectrum  exhibits  nearly  two  hundred  per¬ 
cent  error.  Finally,  we  have  also  computed  the  discrete 
Stokes  spectrum  for  a  wide  variety  of  mesh  resolutions 
and  found  similar  results  to  those  reported  here.  This  is 
also  the  case  for  the  other  discretization  schemes  con¬ 
sidered  in  this  paper. 


3.2  The  Non-conforming  Crouzeix-Raviart  Element 

We  proceed  by  conducting  a  discrete  spectrum  analysis 
for  the  non-conforming  Crouzeix-Raviart  element  [11], 
For  this  element,  midside  nodes  are  used  as  DOF  for  the 
discrete  velocity  field.  This  generates  a  non-conforming 
piecewise  linear  approximation.  The  discrete  pressure 
field  is  taken  to  be  constant  on  each  element.  The  DOF 
for  this  element  are  illustrated  in  Fig.  1(b).  As  the  di¬ 
vergence  of  the  discrete  velocity  field  is  piecewise  con¬ 
stant,  this  velocity /pressure  pair  results  in  a  discrete 
velocity  field  which  is  exactly  divergence-free  over  each 
element.  Furthermore,  the  non-conforming  Crouzeix- 
Raviart  element  admits  an  element-wise  mass  conser¬ 
vation  law  in  the  sense  that  the  total  amount  of  mass 
which  flows  out  of  the  element  is  equal  to  the  total 
amount  of  mass  which  flows  in.  However,  the  element 
is  not  strongly  mass-conservative  as  the  discrete  veloc¬ 
ity  field  associated  with  the  element  does  not  exhibit 
normal  continuity  across  element  boundaries.  Like  the 
MINI  element,  the  non-conforming  Crouzeix-Raviart  el¬ 
ement  is  computationally  inexpensive  and  popular  in 


Fig.  3  Stokes  spectrum  for  the  Crouzeix-Raviart  element. 


practice,  though  it  is  generally  less  robust  than  the 
MINI  element.  For  example,  the  discretized  momentum 
equation  associated  with  the  non-conforming  Crouzeix- 
Raviart  element  is  not  coercive  in  the  presence  of  trac¬ 
tion  boundary  conditions. 

We  have  computed  the  discrete  Stokes  spectrum 
corresponding  to  the  non-conforming  Crouzeix-Raviart 
element  on  a  structured  triangular  mesh  with  2,048  el¬ 
ements  and  plotted  a  comparison  of  the  discrete  spec¬ 
trum  with  the  exact  spectrum  in  Fig.  3.  For  this  calcu¬ 
lation,  N  =  4, 095.  Note  that  all  of  the  discrete  eigen¬ 
values  are  smaller  in  magnitude  than  their  continu¬ 
ous  counterpart.  This  indicates  that,  in  an  unsteady 
flow  simulation,  all  of  the  discrete  modes  will  dissipate 
slower  than  their  continuous  counterparts.  In  a  turbu¬ 
lent  flow  simulation  free  of  turbulence  models,  this  will 
have  the  effect  of  lengthening  the  discrete  inertial  range 
and  pushing  back  the  discrete  viscous  range  to  higher 
wavenumbers.  Furthermore,  while  most  of  the  discrete 
eigenvalues  exhibit  less  than  fifty  percent  error,  less 
than  ten  percent  of  the  discrete  eigenvalues  exhibit  less 
than  ten  percent  error. 

3.3  The  Marker-and-Cell  Scheme 

We  next  consider  the  Marker-and-Cell  (MAC)  finite  dif¬ 
ference  scheme.  The  MAC  scheme  is  one  of  the  old¬ 
est  schemes  for  incompressible  fluid  flow  and  was  origi¬ 
nally  introduced  by  Harlow  and  Welch  to  simulate  free- 
surface  problems  [17].  In  the  MAC  scheme,  the  veloc¬ 
ity  field  is  defined  on  cell  faces,  the  pressure  field  is 
defined  on  cell  interiors,  and  finite  differences  are  uti¬ 
lized  to  discretize  the  momentum  equation  on  cell  faces 
and  the  continuity  equation  on  cell  interiors.  The  DOF 
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for  this  scheme  are  illustrated  in  Fig.  1(c).  The  MAC 
scheme  has  been  shown  to  be  equivalent  to  a  finite  el¬ 
ement  scheme  [24],  a  mixed  finite  element  method  in 
vorticity-velocity-pressure  formulation  with  special  nu¬ 
merical  quadrature  [16],  and  a  particular  divergence- 
conforming  discontinuous  Galerkin  method  [22],  The 
MAC  scheme  satisfies  a  mass  conservation  statement  on 
each  cell  by  construction,  and  this  conservation  prop¬ 
erty  has  resulted  in  enhanced  robustness  for  problems 
where  mass  conservation  is  of  critical  importance  such 
as  coupled  fluid-transport  [23]. 

We  have  computed  the  discrete  Stokes  spectrum 
corresponding  to  the  MAC  scheme  on  a  structured  rect¬ 
angular  mesh  with  4,096  elements  and  plotted  a  com¬ 
parison  of  the  discrete  spectrum  with  the  exact  spec¬ 
trum  in  Fig.  4.  For  this  calculation,  N  =  4,095.  Note 
that  all  of  the  discrete  eigenvalues  are  smaller  than  their 
continuous  counterpart  as  was  the  case  for  the  non- 
conforming  Crouzeix-Raviart  element.  Hence,  they  will 
suffer  from  the  same  issues  as  the  Crouzeix-Raviart  el¬ 
ement.  Namely,  in  a  turbulent  flow  simulation  without 
the  introduction  of  turbulence  models,  the  discrete  in¬ 
ertial  range  will  be  artificially  lengthened.  However,  for 
the  most  part,  the  discrete  spectrum  corresponding  to 
the  MAC  scheme  is  better  resolved  than  that  of  the 
Crouzeix-Raviart  element.  Approximately  ten  percent 
of  the  discrete  eigenvalues  exhibit  less  than  ten  percent 
error. 


3.4  The  Taylor-Hood  Family 

We  continue  by  conducting  discrete  spectrum  analyses 
for  the  Taylor-Hood  family  of  Stokes  elements  [18].  For 
a  given  polynomial  degree  k  >  2,  the  velocity  held  is 


(b) 


Fig.  5  Degrees  of  freedom:  (a)  the  second-order  Taylor- 
Hood  pair,  (b)  the  Q2/Pi  pair,  (c)  the  second-order  Raviart- 
Thomas  pair.  Velocity  DOF  are  on  the  left  and  pressure  DOF 
are  on  the  right. 


approximated  by  continuous  piecewise  tensor-product 
polynomials  of  degree  k  and  the  pressure  held  is  approx¬ 
imated  by  continuous  piecewise  tensor-product  polyno¬ 
mials  of  degree  k  —  1.  We  denote  this  velocity/pressure 
pair  as  Qk/Qk-i ■  The  DOF  for  the  second-order  Taylor- 
Hood  pair  Q2/Q1  are  illustrated  in  Fig.  5(a).  The  second- 
order  Taylor-Hood  pair  is  one  of  the  most  popular  fi¬ 
nite  elements  for  incompressible  how  simulation.  It  is 
balanced  in  the  sense  that  the  discrete  velocity  held  is 
second-order  accurate  in  the  Hi-norm  and  the  discrete 
pressure  held  is  second-order  accurate  in  the  L2-norm. 
Mixed  finite  elements  which  are  unbalanced  typically 
suffer  from  suboptimal  convergence. 

For  k  =  2,3,4,  5,  we  have  computed  the  discrete 
Stokes  spectrum  corresponding  to  the  Qk/Qk-i  ele¬ 
ment  on  a  structured  rectangular  mesh  and  plotted 
a  comparison  of  the  discrete  spectrum  with  the  ex- 


6 


John  A.  Evans,  Thomas  J.R.  Hughes 


Fig.  6  Stokes  spectrum  for  the  Taylor-Hood  family. 


act  spectrum  in  Fig.  6.  For  each  polynomial  degree,  we 
chose  a  sufficiently  fine  mesh  such  that  N  ss  4000.  Note 
that  the  initial  portion  of  the  discrete  spectrum  im¬ 
proves  dramatically  with  increased  polynomial  degree. 
This  is  expected  from  classical  a  priori  error  estimates. 
Asymptotically,  discrete  eigenvalues  converge  like  h2k 
for  the  Taylor-Hood  family.  However,  the  entire  dis¬ 
crete  spectrum  does  not  improve  with  increased  poly¬ 
nomial  degree.  In  fact,  the  upper  part  of  the  discrete 
spectrum  worsens.  This  indicates  that  the  high  modes 
of  the  discrete  Stokes  system  have  no  approximability 
and  these  modes  get  worse  with  increased  polynomial 
degree.  This  is  analogous  to  what  has  been  observed  for 
hp- finite  elements  in  the  context  of  structural  dynamics 
[21].  Furthermore,  it  should  be  noted  that  many  of  the 
discrete  eigenvalues  corresponding  to  the  Q2/Qi  pair 
are  smaller  than  their  continuous  counterpart. 

3.5  The  Qk/Pk-i  Family 

We  next  consider  the  Qk/Pk-i  family  of  Stokes  ele¬ 
ments  [3].  In  this  discretization  scheme,  the  velocity 
field  is  approximated  using  continuous  piecewise  tensor- 
product  polynomials  of  degree  k  and  the  pressure  field  is 
approximated  using  discontinuous  piecewise  polynomi¬ 
als  of  degree  k—  1.  The  DOF  for  the  second-order  mem¬ 
ber  of  the  family  are  illustrated  in  Fig.  5(b).  Like  the 
Taylor-Hood  pairs,  the  Qk/Pk-i  pair  is  balanced.  Fur¬ 
thermore,  as  the  discrete  pressure  field  contains  piece- 
wise  constants,  the  discrete  velocity  field  arising  from 
a  Qfc/Tfc-i  discretization  is  mass-conservative  element- 
by-element.  Like  the  MAC  scheme,  this  conservation 
property  has  resulted  in  enhanced  robustness  for  prob¬ 
lems  where  mass  conservation  is  of  importance. 


For  k  =  2,3,4,  5,  we  have  computed  the  discrete 
Stokes  spectrum  corresponding  to  the  Qk/Pk-i  on  a 
structured  rectangular  mesh  and  plotted  a  comparison 
of  the  discrete  spectrum  with  the  exact  spectrum  in  Fig. 
7.  For  each  polynomial  degree,  we  chose  a  sufficiently 
fine  mesh  such  that  N  ss  4000.  Note  immediately  that 
the  discrete  spectra  associated  with  these  discretiza¬ 
tions  are  somewhat  erratic  and  unpredictable.  Like  the 
Taylor-Hood  family  of  discretizations,  the  initial  por¬ 
tion  of  the  discrete  spectrum  improves  with  increasing 
polynomial  degree  and  the  upper  part  worsens.  More¬ 
over,  as  the  polynomial  degree  is  increased,  the  mid¬ 
dle  portion  of  the  discrete  spectrum  is  driven  to  be 
smaller  and  smaller.  This  being  said,  for  a  given  poly¬ 
nomial  degree  k,  the  discrete  spectrum  associated  with 
the  Qfc/Pfc_i  finite  element  pair  is  generally  more  ac¬ 
curate  than  the  discrete  spectrum  associated  with  the 
Qk/Qk-i  pair. 

3.6  Divergence-conforming  DG  Elements 

We  proceed  by  conducting  discrete  spectrum  analy¬ 
ses  for  divergence-conforming  discontinuous  Galerkin 
(DG)  schemes  which  arise  when  one  employs  a  Raviart- 
Thomas  velocity/pressure  pair  in  conjunction  with  the 
symmetric  interior  penalty  Galerkin  method  to  enforce 
tangential  continuity  across  element  faces  [8,9,22],  Such 
schemes  have  become  increasingly  popular  in  recent 
years  as  they  result  in  discrete  velocity  fields  which 
are  pointwise  divergence-free.  In  the  context  of  Navier- 
Stokes  flow,  exact  satisfaction  of  the  incompressibility 
constraint  guarantees  that  divergence-conforming  DG 
schemes  automatically  conserve  momentum  and  are  en¬ 
ergy  stable.  The  DOF  associated  with  the  second-order 
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Fig.  8  Stokes  spectrum  for  divergence-conforming  DG  ele¬ 
ments.  Top:  global  view,  bottom:  zoomed  in  view. 


Raviart-Thomas  pair  are  depicted  in  Fig.  5(c).  It  should 
be  noted  that  in  the  context  of  Stokes  flow,  the  classical 
second-order  Raviart-Thomas  pair  is  only  first-order  ac¬ 
curate  as  the  discrete  space  of  velocity  fields  associated 
with  the  second-order  Raviart-Thomas  pair  is  complete 
only  up  to  linear  polynomials.  For  this  reason,  we  will 
designate  the  divergence-conforming  DG  scheme  em¬ 
ploying  a  kth- order  Raviart-Thomas  velocity/pressure 
pair  as  ( k')th  =  (k  —  l)*/t-order. 

For  k'  =  1,2, 3, 4, 5,  we  have  computed  the  dis¬ 
crete  Stokes  spectrum  corresponding  to  the  (fc'^-order 
divergence-conforming  DG  scheme  on  a  structured  rect¬ 
angular  mesh  and  plotted  a  comparison  of  the  discrete 
spectrum  with  the  exact  spectrum  in  Fig.  8.  For  each 
polynomial  degree,  we  chose  a  sufficiently  fine  mesh 
such  that  N  «  4000,  and  we  selected  the  interior  penalty 
parameter  as  4(fc'  +  1)2.  Note  from  the  figure  that  while 
the  initial  portion  of  the  discrete  spectrum  improves 


Fig.  9  Degrees  of  freedom:  divergence-conforming  B-splines. 
Velocity  DOF  are  designated  by  arrows  and  pressure  DOF 
are  designated  by  points. 

with  increased  polynomial  degree,  the  upper  part  of 
the  spectrum  catastrophically  worsens.  Furthermore, 
we  found  the  quality  of  the  discrete  spectrum  degraded 
with  increasing  interior  penalty  parameter.  This  indi¬ 
cates  that  one  should  choose  the  interior  penalty  pa¬ 
rameter  intelligently  not  only  for  the  purposes  of  nu¬ 
merical  stability  but  also  for  numerical  accuracy.  Un¬ 
fortunately,  this  fact  is  not  often  elaborated  in  the  lit¬ 
erature. 

3.7  Divergence-conforming  B-splines 

We  finish  by  conducting  discrete  spectrum  analyses  for 
divergence-conforming  B-spline  approximations  of  the 
Stokes  eigenproblem  [5,12].  These  approximations  are 
motivated  by  the  recent  theory  of  isogeometric  discrete 
differential  forms  [6, 7]  and  may  be  interpreted  as  smooth 
generalizations  of  Raviart-Thomas  elements.  As  these 
approximations  are  smooth,  they  can  be  directly  uti¬ 
lized  in  the  Galerkin  solution  of  viscous  flows  without 
resorting  to  a  DG  technology  in  contrast  with  classi¬ 
cal  Raviart-Thomas  elements.  Recent  work  has  shown 
these  approximations  are  much  more  accurate  and  sta¬ 
ble  than  classical  methods  for  Stokes  and  Navier-Stokes 
flows  [13-15].  Furthermore,  in  the  context  of  Navier- 
Stokes  flow,  these  approximations  admit  balance  laws 
for  momentum,  energy,  vorticity,  enstrophy,  and  he- 
licity  as  they  satisfy  the  incompressibility  constraint 
exactly  [15].  Since  divergence-conforming  B-spline  ap¬ 
proximations  exhibit  additional  smoothness  in  compar¬ 
ison  with  classical  finite  elements,  it  is  not  possible  to 
attribute  DOF  to  element  quantities.  Instead,  DOF  can 
be  identified  to  cells  and  edges  on  a  global  data  struc¬ 
ture  called  a  control  mesh.  For  reference,  see  Chapter  4 
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Fig.  10  Stokes  spectrum  for  divergence-conforming  13- 
splines. 


of  [12]  or  Chapter  2  of  [10].  We  have  depicted  a  standard 
DOF  pattern  on  a  patch  in  Fig.  9.  Note  the  similarity 
of  the  DOF  structure  with  that  of  the  MAC  scheme. 

For  k'  =  1,2, 3, 4, 5,  we  have  computed  the  dis¬ 
crete  Stokes  spectrum  corresponding  to  the  (fc'^-order 
divergence-conforming  B-spline  scheme  on  a  structured 
rectangular  mesh  with  4, 096  elements  and  plotted  a 
comparison  of  the  discrete  spectrum  with  the  exact 
spectrum  in  Fig.  10.  For  these  calculations,  N  =  4, 095. 
Note  that  the  discrete  spectrum  is  quite  well-resolved 
for  all  polynomial  degrees  and  that  the  entire  discrete 
spectrum  improves  markedly  with  increasing  polyno¬ 
mial  degree.  The  discrete  spectrum  corresponding  to 
the  first-order  B-spline  discretization  is  remarkably  more 
accurate  than  the  discrete  spectra  corresponding  to  the 
other  first-order  methods  in  this  paper  (the  MINI  ele¬ 
ment,  the  Crouzeix-Raviart  element,  the  MAC  scheme, 
and  the  first-order  divergence-conforming  DG  scheme), 
and  the  discrete  spectrum  corresponding  to  the  second- 
order  B-spline  discretization  is  dramatically  more  ac¬ 
curate  than  the  discrete  spectrum  associated  with  any 
other  method  considered  in  this  paper.  In  fact,  ap¬ 
proximately  eighty  five  percent  of  the  spectrum  ex¬ 
hibits  less  than  ten  percent  error  for  the  second-order 
B-spline  discretization,  and  approximately  eighty  five 
percent  of  the  spectrum  exhibits  less  than  four  per¬ 
cent  error  for  the  fifth-order  B-spline  discretization.  We 
should  finally  mention  that  all  of  the  discrete  eigenval¬ 
ues  associated  with  a  divergence-conforming  B-spline 
approximation  necessarily  satisfy  >  Xn  since  the 
approximation  scheme  is  a  conforming  interior  approx¬ 
imation  technique  (i.e. ,  divergence-conforming  B-spline 
schemes  deliver  discrete  velocity  fields  which  are  point- 
wise  divergence- free). 


4  Conclusions 

In  this  paper,  discrete  spectrum  analyses  were  conducted 
for  a  selection  of  mixed  discretization  schemes  for  the 
Stokes  eigenproblem.  These  analyses  revealed  that  clas¬ 
sical  schemes  such  as  the  MINI  element,  the  non-conforming 
Crouzieux-Raviart  element,  and  the  Marker-in-Cell  scheme 
exhibit  poorly  resolved  discrete  spectra,  and  they  fur¬ 
ther  revealed  that  the  upper  part  of  the  discrete  spec¬ 
trum  degraded  under  degree  elevation  for  Taylor-Hood 
elements,  Qfc/Pfc_i  elements,  and  divergence-conforming 
DG  elements.  This  is  concerning  as  the  eigenvalues  of 
the  Stokes  eigenproblem  correspond  to  the  dissipation 
rates  of  the  natural  modes  of  unsteady  Stokes  flow  as 
well  as  the  squares  of  the  natural  frequencies  of  linear 
incompressible  elastodynamics.  Consequently,  a  scheme 
with  favorable  resolution  properties  should  be  charac¬ 
terized  by  a  discrete  spectrum  which  closely  mimics  the 
exact  spectrum.  In  contrast  with  the  other  methods 
studied  in  this  paper,  it  was  discovered  that  the  dis¬ 
crete  spectrum  associated  with  divergence-conforming 
B-spline  approximations  is  quite  well-resolved  for  all 
polynomial  degrees  and  that  the  entire  discrete  spec¬ 
trum  improves  markedly  with  increasing  polynomial  de¬ 
gree.  This  motivates  further  study  of  divergence-conforming 
B-spline  approximations  in  the  context  of  fluid  flow  sim¬ 
ulation  and  incompressible  elastodynamics. 
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